TK!2
=v
#1
:n

:s
b
:c
Newton's Method for Roots of Functions
#3
:n
fun
:v
'poly4
:s
i
:c
name of function
#4
:n
dfun
:v
'dpoly4
:s
i
:c
name of symbolic derivative
#5
:n
x0
:v
.1
:s
i
:c
guess value
#6
:n
eps
:v
.000001
:s
i
:c
tolerance 
#7
:n
xS
:v
2.000000004719277
:s
o
:c
solution returned by NewtonS
#8
:n
xN
:v
2.000000023234354
:s
o
:c
solution returned by NewtonN
#9
:n
dx
:v
.0001
:s
i
:c
increment for numer. differentiation
#10
:n
y
:s
b
:c
) used to generate graph of f(x)
:a
y
#11
:n
x
:v
0
:s
i
:c
)
:a
x
=u
=l
#1
:n
x
:v
-2.375
-2.25
-2.125
-2
-1.875
-1.75
-1.625
-1.5
-1.375
-1.25
-1.125
-1
-.875
-.75
-.625
-.5
-.375
-.25
-.125
0
.125
.25
.375
.5
.625
.75
.875
1
1.125
1.25
1.375
1.5
1.625
1.75
1.875
2
2.125
2.25
2.375
E
#2
:n
y
:v
7.614
4.316
1.813
0
-1.219
-1.934
-2.23
-2.188
-1.879
-1.371
-.726
0
.758
1.503
2.199
2.812
3.317
3.691
3.922
4
3.922
3.691
3.317
2.812
2.199
1.504
.758
0
-.726
-1.371
-1.879
-2.188
-2.23
-1.934
-1.219
0
1.813
4.316
7.614
E
=f
#1
:n
poly4
:c
function example
:t
Rule
:a
x
:v
y
#2
:n
dpoly4
:c
dy/dx by symbolic differentiation
:t
Rule
:a
x
:v
dy%dx
#3
:n
NewtonS
:c
Newton's method, symbolic derivatives
:t
Procedure
:i
FUN,dFUN,x0,eps
:o
x
#4
:n
NewtonN
:c
Newton's method, numeric differentiation
:t
Procedure
:i
FUN,x0,dx,eps
:o
x
#1
:n
poly4
:z
y = poly(x,1,0,-5,0,4)   " i.e.  y = x^4 - 5*x^2 + 4
/E
E
#2
:n
dpoly4
:z
dy%dx = poly(x,4,0,-10,0)   " i.e.  dy%dx = 4*x^3 - 10*x
/E
E
#3
:n
NewtonS
:z
x:= x0
Loop:
   y:= apply(FUN,x)
   if abs(y) < eps  return
   x:= x - y/apply(dFUN,x)
goto Loop
/E
E
#4
:n
NewtonN
:z
x:= x0
Loop:
   y:= apply(FUN,x)
   if abs(y) < eps  return
   x:= x - y/((apply(FUN,x*(1+dx)) - y) / (x*dx))
goto Loop
/E
E
=r
#1
:r
"        >>>>  Newton's Method for Finding Roots of Functions  <<<<
:s
C
#3
:r
xS = NewtonS(fun,dfun,x0,eps)      "   cancel when generating graph of y=f(x)
#4
:r
xN = NewtonN(fun,x0,dx,eps)        "     "     "       "        "   "    "
#5
:r
y = apply(fun,x)                   " uncancel  "       "        "   "    "
:s
C
#7
:r
" Newton's (also called Newton-Raphson) method is the basis of TK's Iterative
:s
C
#8
:r
" Solver, which can easily handle one or more unknown variables.  However,
:s
C
#9
:r
" occasionally 
:s
C
#10
:r
"    you may not want to use the Iterative Solver (or you need it in
:s
C
#11
:r
"     another part of your model),
:s
C
#12
:r
"    an unknown you want to evaluate cannot be solved directly,
:s
C
#13
:r
"    the function is well behaved (i.e. it is monotonic in the interval
:s
C
#14
:r
"     of your interest) and a reasonable guess is available or can be
:s
C
#15
:r
"     calculated from other variables.
:s
C
#16
:r
" If this is so, you can build one or more Newton's root finders within your
:s
C
#17
:r
" model.  This technique is used in the inverse Chi-Square and F cumulative
:s
C
#18
:r
" distribution functions in the model STAT\DISTRIB and in the function RATEL
:s
C
#19
:r
" in the model FINANCE\FILCAL2.  It is discussed at length in the model
:s
C
#20
:r
" DIFFINT\SYMDIFF.
:s
C
#22
:r
" This model presents another example in two versions - with symbolic
:s
C
#23
:r
" derivatives (which is preferred when the derivatives are available) and
:s
C
#24
:r
" with numeric differentiation.  The sample function is a 4th order poly-
:s
C
#25
:r
" nomial  y = x^4 - 5*x^2 + 4  with roots  -2, -1, 1, 2 (see plot "poly4").
:s
C
#27
:r
" You can experiment with these function and see that they work fast and
:s
C
#28
:r
" well, although the resulting values depend on the initial guess.  If you
:s
C
#29
:r
" try guessing  x0=0  the solution will be terminated by a "Division by zero"
:s
C
#30
:r
" condition (reason: the tangent to the poly4 curve at that point is
:s
C
#31
:r
" horizontal -- a situation which has to be avoided when using Newton's
:s
C
#32
:r
" method).
:s
C
=p
#1
:n
poly4
:t
Line chart
:d
EGA
:i
"y  =  x^4  -  5 * x^2  +  4"
:s
Yes
:z
Both
:g
No
:l
None
:h
"x"
:v
"y"
:x
x
:u
#1
:y
y
:f
1
:p
Lines
:c
"*"
E
=t
=n
=g
#1
:c
.000001
:i
10
:v
Yes
:l
60
:f
66
:m
0
:w
80
:b
No
:a
Yes
:r
Yes
:h
No
:t
1
:n
No
:d
LPT1
%Tv,f,0,5,0,1,0
%Or,f,0,1,0,1,0
